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A.  Calculation  of  Wave  Resistance  and  Sinkage  by  Rankine-Source  Method 

1.  Introduction 

The  Rankine  source  was  used  first  by  Gadd  [1]  for  the  free-surface 
flow  problem.  Later  Dawson  derived  a  more  elegant  formula  and  presented  many 
numerical  results  which  may  be  considered  a  milestone  in  ship  hydrodynamics 
[2,3]. 

One  of  the  significant  advantages  of  the  method  is  its  generality 
and  the  simple  form  of  the  Green  function.  When  the  double-hull  linearized 
free-surface  condition  is  used,  it  can  be  extended  to  the  viscous  flow  problem 
(Mori  X  Nishimoto  [4]).  Not  only  the  wave  resistance  but  also  the  moment  and 
the  force  acting  on  a  ship  hull  can  be  calculated.  It  enables  us  to  attain 
our  final  goal  to  calculate  the  wave  resistance  of  an  unrestrained  ship  form, 
including  viscous  effects. 

In  the  present  interim  report,  calculations  are  limited  to  the  case 
of  potential  flow  under  the  restrained  condition.  This  is  because  there  are 
still  several  points  which  should  be  clarified  in  the  numerical  techniques. 
The  inclusion  of  viscosity  and  the  iterative  calculation  to  realize  the  trim- 
and  sinkage-free  condition  are  left  for  the  second  half  of  the  present  work. 

2.  Basic  equations 

Let  <j>(x,y,z)  be  the  perturbation  velocity  potential  at  P(x,y,z)(see 
Fig.  1  for  definitions).  By  Green's  theorem,  4  can  be  expressed  in  the  form 

4-*(P )  =  -JJ  o H (“  +  "pr ) dS  -  }J  oF”dS  +  4n^,  (1) 

SH  SF 

where  SH  denotes  the  hull  surface,  Sp  the  free  surface,  a  source 
distribution  over  S^,  and  a  source  distribution  over  Sp,  (see  Appendix). 


r2=  (x  -  x')2  +(y  -  y')2  +  (z-z1)2 
r'2=  (x  -  x')2  +  (y  -  y')2  +  (z  +  z’)2, 


(2) 


.1 
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(x,y,z)  are  the  coordinates  of  a  fixed  point  P,  (x'.y'.z')  those  of  the 
integration  points.  41  is  the  contribution  from  the  integration  over  the  far- 

00 

field  surface  S0,  given  by 


*♦.- 11  <♦  h  <r *  &  -  (f  *  & «•  <3) 

J3 

In  the  domain  of  VQ,  which  is  partially  surrounded  by  S0  (see  Fig. 
1),  we  can  assume  that  the  velocity  potential  4  satisfies  the  linearized  free- 
surface  condition 


K 


Ml 

O  3  Z 


(4) 


2 

where  KQ  =  U  /g,  g  is  the  gravity  acceleration,  and  U  is  the  ship  speed.  By 
making  use  of  the  Havelock  Green  function  G  which  satisfies  the  linearized 
free-surface  condition,  application  of  Green's  formula  gives 


-  1  <*'!v-G&>  *'  *  11  <♦'& 
0 

0  0 


G  &  >  *• 


(5) 


where  L  is  the  intersection  of  S  and  the  free  surface. 
0  0 


Because  and  34.' /a x  should  be  continuous  on  S0,  the  addition  of  Eq.  (5)  to 
Eq.  (3)  yields 


4tt4>  = 


4j. 

0  L. 


(♦ 


3  G 
3  x' 


34> 1 


3  X 


■->  <t y-ii  (* '  - H  >  ds> 


3  4> 1 


(6) 


where 


H  =  G  +  (-p  +  “f ) . 


(7) 


Due  to  the  exponential  decay  of  H  in  the  depthwise  directiion,  the  integration 
over  S0  may  be  limited  close  to  the  free  surface. 

Eq.  (1)  has  three  unknowns;  oH,  op  and$'.  They  can  be  determined  so 
that  +  satisfies  the  free-surface  condition,  the  hull-surface  condition  and  an 
additional  condition  on  SQ.  This  is  accomplished  iteratively.  We  write  the 
velocity  in  the  form 
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q  =  V  7V  (8) 

where 

qo=  iU  +  7*o  +  V  {9) 

and 

♦o'  '  t  1L  °H<r  *  7T>ds  <10> 

*w=  -  V  ii  °F  r  dS  +  *»  (U) 

bF 

Here  1  is  the  unit  vector  in  the  x-direction,  and  qv  a  viscous  velocity 
vector. 

At  first,  we  determine  oH  in  order  to  satisfy  the  hull-surface 
condition  for  a  given  4>w.  At  the  first  iteration,  Op  and  ^  are  assumed 
zero.  Thus  q0,  at  the  first  iteration,  gives  the  double-hull  flow  itself. 
Then  op  is  determined  so  as  to  satisfy  the  free-surface  condition  and  a  proper 
downstream  condition  imposed  on  S0* 

The  free-surface  condition,  which  is  linearized  based  on  the  flow 
field  of  q0,  is 


^o^o^wtV*  qo<Wwx+  qoy^wy)  +  ^wz  "  qoqot  * 

on  z  -  o 

where  qQ  =  |qQ|  Z=Q  and  subscripts  l ,  x  and  y  imply  differentiation  with 
respect  to  the  indicated  variables.  Here  t  denotes  arclength  along  a 
streamline  on  z  =  o. 

The  condition  imposed  on  Sq  is  that  <j>  ’  must  be  matched  with  a 
solution  of  the  Laplace  equation  which  satisfies  the  linearized  free-surface 
condition  and  the  radiation  condition.  In  the  present  calculation, 
however,  ^  is  neglected  by  choosing  the  computing  domain  and  the  finite- 
difference  scheme  properly  according  to  the  pilot  computations  given  in  [4]. 

The  transformation  of  Eq.  (12)  into  a  finite-difference  equation 
provides  a  set  of  simultaneous  equations  for  the  unknown  variable  Op.  Its 
precise  expression  can  be  found  in  [4]. 
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Once  the  velocity  potential  is  determined,  the  pressure  on  the  hull 
is  given  by 


P  =  7  p  (U2  -  |q|2)  -  pg6H,  (13) 

where  q  is  the  total  velocity  vector,  p  is  the  density  of  fluid,  and  6H  the 
head  loss  whose  gradient  is  given  by 


c^SH  =  m  x  q  -  vV^q,  (14) 

where  «  is  the  vorticity  vector  and  v  is  the  kinematic  viscosity.  Then,  the 
pressure  resistance,  Rp,  the  sinkage  force,  Rs,  and  the  trim-by-stern  moment, 
Mt,  are  given  as  follows: 


Rp  ’-11  SH  pnxdS' 


Rs  =  JJ  pnzdS, 


(15) 


Mt  =  iisH  p{nz(x  ‘  xo}  ‘  nx(z  "  zo)}dS’ 

where  (x0,  zQ)  are  the  coordinates  of  the  center  of  buoyancy,  and  (nx,  nz)  are 
the  x-  and  z-components  of  the  unit  outward  normal  vector  on  SH. 


3.  Numerical  calculation  and  discussions 

Two  ships  are  chosen  for  the  present  calculation;  an  Inuid  model  M- 
21  [5]  and  the  Wigley  parabolic  model  whose  principal  dimensions  are  shown  in 
Table  1.  The  calculations  are  limited  to  the  inviscid  case  and  the  first 
iteration  is  carried  out.  Because  M-21  is  a  model  generated  by  the  streamline 
tracing  method,  is  exactly  known,  thus  avoiding  some  numerical  errors. 

Fig.  2  shows  the  discretization  of  the  free  surface  for  M-21.  For 
the  present  calculation,  the  computing  domain  is  chosen  as  -  1.5  <  2x/L  <  5.0 
and  0  <  2x/l  <  4.5,  where  L  is  the  ship  length. 

The  first  derivatives  appearing  on  the  l.h.s.  of  Eq.  (12)  are 
approximated  by 


Fi 


in  <-  Fi-3  +  6Fi. 


2-15Fi_1+  10F i ) , 


■  3r.  7V*i> 


(16) 


t 
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where  h  is  the  grid  size,  while  those  on  the  r.h.s.  (knowns)  are  obtained  by  a 
simple  centered  difference.  Eq.  (17)  is  the  form  which  keeps  terms  up  to  the 
fourth  derivatives  in  the  Taylor  expansion  but  drops  the  third.  It  has  been 
determined  in  [4]  that  the  finite  difference  of  Eq.  (17)  is  better  than  other 
forms;  it  avoids  downstream  reflection  and  dampens  waves,  making  It  possible 
to  satisfy  the  radiation  condition. 

The  upstream  condition  is  satisfied  by  requiring  that  the 
perturbation  velocity  be  zero  on  the  three  upstream  columns  of  panels 
(shaded  in  Fig.  2). 

Fig.  3  shows  the  resulting  free-surface  source  distribution  c’ong 
the  2nd  and  3rd  rows  of  panels  (j  =  2,3).  For  comparison,  the  results 
obtained  by  changing  the  computing  domain  from  -1.5  <  2x/L  <  5.0  to  -1.5  < 
2x/L  <  3.0,  are  also  shown.  The  results  are  quite  similar.  This  does  not 
always  guarantee  the  exact  satisfaction  of  the  radiation  condition,  but  it  can 
be  safely  concluded  that  the  truncation  of  the  computing  domain  does  not 

affect  the  results  greatly  if  the  selected  range  extends  at  least  two  ship 

lengths  upstream  and  downstream  from  the  ship  form  and  the  finite  difference 
procedure  is  properly  chosen. 

Fig.  4a  shows  the  comparisons  of  perturbation  velocity  components, 
u,v  and  w,  in  the  x-,  y-  and  z-directions  respectively.  The  measurement  is 
carried  out  at  2z/L  =  -0.02  (just  beneat  the  free  surface)  while  the  computed 
results  are  on  z=o.  The  values  of  u  are  in  good  agreement  except  near  the  bow 
and  stern.  The  discrepancies  near  the  bow  may  be  due  to  the  approximation 

that  <t>  does  not  satisfy  the  hull-surface  condition,  while  those  near  the 

stern  may  be  attributable  to  the  effects  of  viscosity.  (It  is  intended  to 
include  the  viscous  effects  in  subsequent  calculations).  The  comparisons  of 
the  v-  and  w-components  indicates  the  need  for  a  slight  refinement  of  the 
numerical  techniques.  As  seen  in  Fig.  4b,  however,  the  comparisons  at  a 
deeper  position,  2z/L  =  -  0.1,  shows  much  better  agreement. 

In  Fig.  5  the  calculated  wave  resistance  and  sinkage  are  compared 
with  measured  results.  The  measurement  of  resistance  was  carried  out  under 
the  sinkage-free  condition  (trim-fixed).  The  wave  resistance  is  obtained  by 
subtracting  the  viscous  resistance  from  the  total  resistance.  The  sinkage, 
2s/L,  was  calculated  by  means  of  the  static  relation,  given  by 
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The  calculated  wave  resistance  is  smaller  than  that  measured.  This 
discrepancy  may  come  mainly  from  the  inadequate  agreement  in  velocity 
(pressure)  observed  near  the  bow,  especially  in  the  u-component.  Another 
possible  reason  is  the  neglect  of  sinkage  in  the  calculation.  The  neglect  of 
the  viscosity  is  also  a  possible  reason,  for  M-21  is  an  optimized  hull  and  the 
inclusion  of  viscosity  may  increase  the  wave  resistance. 

The  disagreement  in  the  sinkage  may  be  due  to  the  same  reasons.  It 
should  be  remembered  that  the  calculated  result  shown  here  is  that  of  only  the 
first  iteration.  Therefore,  we  can  expect  much  better  agreement  in  successive 
iterations  where  the  sinkage  of  ship  is  taken  into  account. 

In  the  case  of  the  Wig  ley  model,  the  hull  surface  is  divided  into  24 
(lengthwise)  x  5  (draftwise)  panels.  The  computing  domain  is  chosen  as  -1.5  < 
2x/L  <  3.0  and  0  <  2y/L  0.45,  and  divided  into  318  panels.  The  calculated 
results  of  wave  resistance  and  sinkage  are  compared  in  Fig.  6.  Experimental 
data  are  those  of  a  4.0m  model  measured  at  ISR  [6],  The  wave  resistance  was 
measured  with  the  model  restricted  in  sinkage,  but  free  to  trim.  Dawson's 
results  [3]  are  also  shown  in  the  figure. 

The  discrepancy,  observed  in  the  wave  resistance  curve  ac  the  higher 
speed  range,  may  be  due  to  the  use  of  the  double-hull  linearized  free-surface 
condition.  The  difference  between  the  present  results  and  Dawsons 's  may  come 
mainly  from  the  numerical  method;  he  solved  Eq.  (1)  directly  to  satisfy  both 
the  hull-surface  condition  and  the  free-surface  condition  simultaneously. 

The  predicted  sinkage  is  smaller  in  this  case  also.  This  is 
probably  die  to  the  same  causes  as  that  for  M-21. 

Several  kinds  of  discretization  of  the  hull  surface  have  been 
compared  for  the  integration  of  pressure  over  the  hull.  These  do  not  make  any 
significant  changes  if  the  source  distributions,  oc  and  au,  are  unchanged.  A 
possible  improvement  may  be  obtained  by  making  use  of  the  Lagally  theorem. 


4.  Concluding  remarks 

The  Rankine  source  method  is  applied  to  the  calculations  of  the  wave 
resistance  and  the  sinkage  of  the  Inuid  model  M-21  and  the  Wigley  model.  Some 
improvements  in  numerical  schemes  and  techniques  may  be  still  necessary,  e.g. 
the  inclusion  of  the  far  downstream  contribution,  to  carry  out  a  successive 
iteration  and  so  on.  It  can  be  safely  concluded,  however,  that  the  method 
works  well  and  is  useful. 

The  author  plans  to  carry  out  extensive  additional  ca leu  lations , 
including  viscous  effects. 
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Appendix  -  Derivation  of  Equation  (1). 

We  consider  three  sub-domains  of  Ve,  Vg  and  which  are  surrounded 
by  surfaces  of  S0,  SQ,  S^,  Sp  S2.  Sg,  S^,  SH  and  Sp  (see  Fig.  1).  P(x,y,z) 
is  a  fix  point  in  Ve. 

By  applying  Green's  formula  to  the  detain  of  Ve,  "Ve  and  Vi,  we  have 


*■*<»-- lls  +S,+S?  (♦IfF-'fII)  dS  +JJsf  (♦  H- K-7&> 

0  1  2  r 


dS 


+  f/c  (4,  |—  -  -  -  -?•*-)  dS,  ( A- 1 ) 
1 J  SH  5  n^  r  r  3nH 


0  -  JJ. 


_  _  (♦ 
so  +  si  + 


ILI  I  24ln 


9  JL  I 


an  r  r  3n^  ^  ^  an  r  r  3nH 


-)  dS 


0  '  '  J1SH  tsj*1  »"H  '  '  *  ,nH)  dS 


♦JJ,  (*hrl-  thr*)*'  <*-V 

jji  n  n 


(A-3 ) 


where  <(>  and  4,.  are  the  velocity  potential  defined  in  Ve  and  Vi  respectively. 
The  definition  of  the  normals  of  n  and  n^  are  shown  in  Fig.  1. 

The  integrands  over  and  "S1  vanish  when  the  upstream  condition  is 
imposed  on  $ .  Then  the  addition  of  Eqs.  (A-l),  (A-2)  and  (A-3)  yields 

<  <♦  *+JJ  C*'h7-}£> 

^o  _ 

s,  ds  +  JJ,  * 

c  2 

+  JJsF  <(♦-♦)  Irr  f  In  (♦  -♦»  « 

*  ^SHt  -  l7br„  ds  *  JJ,  h-f  -  r  rnH  (‘♦-♦i»ds  • 

(A-4) 


It  can  be  expected  that  the  contributions  from  the  integrations  over  S£ 
and  1>2»  the  third  and  fourth  terms  in  Eq.  (A-4),  may  cancel  each  other  when  S£ 
and  S^2  are  taken  deep  enough.  We  assume  that  <t>  is  equal  to  41  at  symmetric 
points;  then  Eq.  (A-4)  can  be  written  as  follows: 
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-JJSc?hr  <♦  -♦>  ^  -  JJs  <♦  h  <7*7')  -  <7 +  (A-5> 

F  0 

(A-5)  can  be  expressed  in  terms  of  source  singularities;  putting 
<j>j,  3/3nH  (4,-^)  =  oH  and3/3n  (<t> )  =  <jp,  we  obtain  Eq.  (1). 
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Table  1  Principal  Paticulars 


M-21 

2.001 

0.2368 

0.1724 

0.0347 

0.6686 


Wigley 

3.000 

0.300 

0.1875 

0.075 

1.329 


Fig.  1  Coordinate  System  and  Definition 


Computing  Domain 


Discretization  of  Free  Surface 


Fig.  4b 


Velocity  Distribution  of  M-21  (2z/L*-0.10) 
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B.  Prediction  of  2-0  Near-Wake  Flow  by  Making  Use  of  Time-Dependent  Vortlclty 
Transport  Equation 

1.  Introduction 

This  study  is  an  extension  of  previous  work  [1]  [2]  where  the  flow  over  a 
ship  stern  is  predicted  by  making  use  of  the  vorticity  transport  equation. 

In  the  first  paper  [1],  an  approximation  is  made  for  the  vorticity 
equation  which  is  basically  the  same  as  the  first-order  boundary- layer 
approximation.  In  the  second  paper  [2],  several  approximations  are  invoked  to 
predict  a  significantly  separated  flow  which  are  then  matched  with  each 
other.  The  full  equation  is  applied  to  a  restricted  region  where  a 
recirculating  flow  is  significant. 

In  both  papers,  the  boundary- layer  equation  is  applied  up  to  the  position 
where  the  flow  is  free  from  separation,  beyond  which  it  is  replaced  by  the 
vorticity  equation.  The  boundary-layer  calculation  provides  initial  or 
boundary  values  for  the  vorticity  equation.  The  total  velocity  field  in  the 
wake  is  given  by  a  sum  of  the  potential  velocity  and  the  induced  velocity  of 
vorticity. 

An  important  difference  from  previous  work  is  that,  in  the  present  work, 
the  vorticity  equation  is  used  without  any  simplifications.  This  means  that 
the  governing  equation  for  the  wake  is  completely  elliptic.  The  process  of 
vorticity  transportation  is  solved  by  the  time-marching  method.  The  k-e  model 
is  invoked  for  the  turbulence  closure  in  the  present  scheme. 

Though  the  vorticity  transport  equation  is  exactly  identical  with  the 
Navier-Stokes  equation,  it  has  several  advantages  over  the  N-S  equation  when 
it  is  applied  to  the  prediction  of  near-wake  flows  of  shiplike  bodies.  First, 
the  velocity  field  can  be  obtained  by  a  simple  sum  of  the  potential  velocity 
and  the  velocity  induced  by  vorticity.  This  removes  the  troublesome 
treatments  of  their  interactions,  for  it  is  not  necessary  to  distinguish  the 
viscous  region  from  the  inviscid.  Another  desirable  feature  is  that  the 
unknown  pressure  term  is  eliminated. 

The  most  important  advantage  is  that  the  computing  domain  of  the 
vorticity  equation  is  confined  to  a  nonzero  vorticity  region.  It  is  more 
definite  and  narrower  than  that  of  the  case  where  the  N-S  equation  is  used. 


The  boundary  values  on  the  downstream  and  the  lateral  terminating  surfaces  are 
definitely  explicit  and  zero.  These  advantages  make  the  numerical  calculation 
much  simpler. 

The  only,  but  important,  difficulty  is  that  of  providing  the  boundary 
values  for  vorticity  on  solid  bodies.  In  the  present  paper,  this  is  treated 
by  a  source  distribution  and  a  vortex  sheet  over  the  body. 

Near  wakes  of  a  flat  plate  and  two  elliptic  cylinders  are  calculated. 


2.  Basic  equations 

Restricting  ourselves  to  2-D  cases,  we  use  Cartesian  co-ordinates  which 
are  attached  to  the  body,  as  shown  in  Fig.  1.  By  the  curl  operation  on  the 
Reynolds  equation,  we  obtain 


3oi  _  f3  (urn) 
3t  ~  "  1  3x 


3 


iMi 

3y 


(IV  +  3U)} 
e'3  x  3y" 


+  2 


.2 

_3 _  ,  ,3jv  3U\i 

3 x3y  ivel3y  "  3x;)  * 
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where  t  is  the  time,  u,v  are  velocity  components  in  x-  ,  y-directions,  u>  the 
vorticity  with  axis  normal  to  the  x-y  plane,  defined  by 


3  v  3  u 

“  “  3X  “  3y 


and  satisfying 


3o>  du>_ 
3  X  3y 


(2) 


(3) 


vg  is  the  equivalent  eddy  kinematic-viscosity  coefficient,  related  to  the 
Reynolds  stress  by 


V  v  +  V 


(4) 


-  =  V. 


.3_V  3_U 

't^3x  +  3y^’ 

where  v  is  the  kinematic  viscosity,  and  u',v 
and  v. 


(5) 


are  turbulence  components  of  u 


The  k-e  equations  used  for  the  turbulence  closure  are 
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The  five  model  constants  are  assigned  the  following  values: 


(6) 


(7) 


(8) 


C0  =  0.09,  Cj 


1.44,  C2=  1.92,  0)<=  1.00,  o£=  1.23. 


Once  the  vorticity  field  is  determined,  the  velocity  at  P(x,y)  is  given  by 


u(x,y) 

v(x.y) 


}  = 


_ _ _ 

(x-x1)  +2(y-y‘ )  2 


r  {y-y  )  t  riv  '  riv  1 
- (x-x 1 )  dx  dy 
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+  {v0(x,y)  > 
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where  Y  denotes  the  nonzero  vorticity  region.  The  derivation  is  given  in  the 
Appendix.  The  first  term  is  the  velocity  induced  by  the  vorticity.  Although 
the  solution  of  Eq.  (1)  and  the  integration  in  Eq.  (9)  are  definitely  confined 
to  Y,  the  induced  velocity  influences  the  flow  beyond  Y.  This  is  the  way  the 
viscid  and  inviscid  interaction  is  taken  into  account  in  the  present  method. 

The  second  term  on  the  r.h.s.  of  Eq.  (9)  is  the  potential  components 
given  by 


u 


3  j> 
0"  3x 


+  u. 


=  i*-  + 

ay 


(10) 


where  4  is  the  total  velocity  potential 
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Uj  and  are  the  velocity  components  due  to  the  image  vorticity  which 
can  be  determined  in  order  for  u  and  v  to  satisfy  the  hull  surface  condition; 
i.e.,  u  =  0,  v  =  0  on  the  hull.  Wu  and  Thompson  [3]  criticized  the  inclusion 
of  these  terms  and  they  used  another  expression  which  does  not  contain  them. 
It  is  true  that  their  expression  is  simpler,  but  it  is  valid  only  when  the 
vorticity  distribution  is  obtained  such  that  the  velocity  satisfies  the 
nor.slip  condition.  As  mentioned  later,  however,  this  is  not  possible  when 
Eqs.  (1)  and  (9)  are  solved  iteratively.  On  the  contrary,  u^  and  v^  can  be 
easily  determined,  as  is  shown  in  the  Appendix,  by  making  use  of  a  source  and 
vortex  distribution.  If  a  constant  distribution  is  assumed  along  a  short 
segment  AB,  Ujand  v1  at  P  are  given  by 


J1 


a,  log  zz  -  ^  <  APB, 
1  AP 


y  i  rp 

v!  =°1  <  APB  -^log 


(11) 


where  Ojand  Yj  are  the  strengths  of  the  source  and  vortex  sheet  respectively 
which  are  determined  so  as  to  satisfy  the  hull-surface  condition. 


3.  Computing  scheme 

We  assume  that,  initially,  a  vorticity  region  over  the  surface  of  a  body 
is  generated  by  a  sudden  motion  of  the  body  to  a  constant  speed  of  U.  The 
flow  is  then  irrotational  except  in  the  vorticity  layer  (boundary  layer). 
Subsequently,  a  wake— nonzero  vorticity  region— is  formed  in  the  vicinity  of 
the  stern  by  the  diffusion  and  convection  of  vorticity.  Thus,  finally,  a 
steady  (or  quasi-steady)  wake  is  formed. 


This  idealization  provides  the  following  computing  scheme  for  this 
Initial  and  boundary-value  problem.  First,  we  carry  out  the  boundary-layer 
calculation  up  to  points  beyond  which  the  flow  is  affected  by  separation. 
Then  we  solve  Eqs.  (4)-(8)  and  Eq.  (1)  with  given  boundary  values  and  a 
prescribed  velocity  field.  Next,  the  velocity  field  is  determined  by  Eqs.  (9) 
and  (10).  The  latter  two  calculations  are  repeated  successively  in  the 
following  order  until  the  solution  converges; 


k  <n+1)  k 

1.  (*}  -  r } 
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+  At 
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2.  >+1>  =  Jn)  +  At(n>fjJn),q<n),v(n+1)),  (13) 

3.  q(^)=  fq(Jn+1),c|n+1),yin+1))  +  7*.  (14) 

where  the  superscript  (n)  refers  to  values  at  t  =  t^.  If  all  the  values 
at  t  =t^  are  known  and  necessary  boundary  values  are  given,  Eqs.  (12)  and 
(13)  can  be  solved  in  the  computing  domain  The  computing  domain 

expands  as  time  elapses  in  order  to  include  the  region  of  nonzero  vorticity 
with  zero  value  on  the  boundary. 

Initial  values,  at  t  =  0,  are  tabulated  in  Table  1.  BgL»  B^,  Bx,  By, 
B$yM  are  boundaries  surrounding  the  computing  domain  V  (see  Fig.  1).  Suffix 
BL  refers  to  values  obtained  by  the  boundary-layer  calculation.  Nonzero 
initial  values  for  k  and  c  are  used  partially  in  V(2x/L  <  1.2).  This  is 
because  the  vorticity  on  does  not  diffuse  so  much  if  k  and  e  are  zero 
everywhere  in  V. 

The  boundary  values  are  shown  in  Table  2. 

As  mentioned  already,  we  cannot  have  any  specified  values  for  vorticity 
on  the  solid  body.  This  is  an  unavoidable  difficulty  of  the  present  method 
where  the  vorticity  equation  is  used.  In  our  scheme  we  supply  this  with 

(n+1)  (n)  ni-\ 

=  qt  /A ,  (15) 

where  A  is  a  short  distance  normal  to  the  hull  surface  and  q^  is  the 
resultant  velocity.  Eq.  (15)  is  obtained  by  making  use  of  the  nonslip 
condition,  but  it  does  not  always  guarantee  the  satisfaction  of  the  nonslip 
condition  at  the  (n+l)-th  step.  0,  and  y,  are  determined  for  the  velocity  at 
the  (n+l)-th  step  to  satisfy  the  hull-surface  condition.  At  the  first  time 
step,  n=0,  the  potential  velocity  is  used;  q^°^=  |v<|>). 

As  mentioned  above,  the  computing  domain  must  be  expanded  indefinitely. 

In  reality,  however,  it  must  be  truncated.  After  Bx  has  arrived  at  some 

position  beyond  which  the  computing  domain  is  not  extended 

2  2 

further,  3  m/3x  =  0  is  imposed  as  a  boundary  condition  on  B  . 
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at,  in  Eqs.  (12)  and  (13),  is  the  time  step.  In  the  present  scheme, 
according  to  Von-Neumann ' s  analysis  [4],  it  is  taken  as 


where  Min  {  }  refers  to  the  minimum  of  the  arguments  defined  in  the  whole 

computing  domain,  ax  and  ay  are  grid  sizes. 

f^  e  and  f  are  finite-difference  expressions  of  the  r.h.s.  of  Eqs.  (1), 
(6)  and  (7).  The  second  upwind-differencing  method  [4]  is  used  for  the  first 

(convection)  terms  on  the  r.h.s.  of  Eqs.  (1),  (6)  and  (7).  On  the  other  hand, 

central  differences  are  used  for  the  diffusion  terms. 

Boundary- layer  calculations  are  carried  out  by  the  integral  method,  with 
Head's  entrainment  formula  used  as  an  auxiliary  equation. 

4.  Numerical  results  and  discussions 
4.1  Flat  plate 

The  main  purpose  of  the  calculation  for  a  flat  plate  is  to  examine  the 
present  scheme  and  numerical  techniques.  All  the  computing  conditions  were 
adjusted  to  those  of  the  experiment  of  Chevray  et  al  [5].  The  overall  length 
L  is  2.4m  and  the  origin  is  shifted  to  the  leading  edge  for  thi;  case.  The 

5 

Reynolds  number  based  on  L  is  6.545  x  10  . 

Results  are  shown  in  Figs.  2-5.  The  positions  where  comparisons  are  made 
with  the  measurements  do  not  exactly  correspond;  the  computed  results  at  x/L  = 
1.025,  1.080,  1.200,  1.600  are  compared  with  experimental  data  at  x/L  =  1.020, 
1.083,  1.208  and  1.625  respectively  (only  the  former  numbers  are  used  later 
on).  The  results  are  those  at  t  =  0.917  L/U,  when  the  computing  domain  has 
expanded  to  x  =  1.72L  and  y  =  0.032L. 

Good  agreement  with  experiment  was  obtained  everywhere  except  far 
downstream  at  x/L  =  1.60.  The  numerical  computations  were  completely  stable. 

Needless  to  say,  because  the  vorticity  beyond  the  terminating  surface  (x 
*  1.72L)  is  not  included,  the  poor  agreement  far  downstream  (x=  1.60L)  can  be 
partially  attributed  to  this  exclusion.  Furthermore,  a  possible  invalidity  of 
the  use  of  the  k-e  model  far  downstream  may  be  mentioned.  Recently,  through 
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precise  comparisons  with  experimental  data,  Ramaprian  et  al  [6]  and  Patel  et 
al  [7]  have  pointed  out  that  the  basic  k-e  model  does  not  lead  to  a 
theoretically  expected  asymptotic  flow  in  a  far  wake  due  to  a  different 
behavior  of  the  intermittency. 

As  far  as  the  near  wake  flow  is  concerned,  however,  we  can  conclude  that 
the  present  scheme  works  well  and  satisfactory  predictions  can  be  obtained. 


4.2  Elliptic  cylinder 

Calculations  were  performed  for  two  elliptic  cylinders  designated  EM-125 
and  EM-200,  having  ratios  b  of  the  major  and  minor  axes  of  0.125  and  0.200 
respectively  [8].  The  velocity  measurements  were  carried  out  by  making  use  of 
a  five-hole  pitot  tube.  An  eight-hole  pipe,  which  has  eight  holes  along  the 
circumference  of  a  circular  pipe,  was  also  used  to  follow  reverse  flows.  The 
referred  Reynolds  numbers  of  EM-125  and  EM-200  are  1.68  x  10^  respectively. 

The  potential  velocity  components  are  given  by 


i$_  -  ii  .  1  +  b  cosh  2s  -  cosh  2n  -  b  sinh  2  £ 

9  9  9  u » 

2(1  +  b  )  cosir  c  -  cos^  n 


_  b (1  +  b) _ sin  2n 

2(1  +  b^)  cosh  c  -  cos^n 


(17) 


where 


x  +  iy  =  /\  +  b^ 


cosh  (£  -  in). 


(18) 


The  boundary- layer  calculations  were  carried  out  by  making  use  of  the 
momentum-thickness  equation  and  the  entrainment  equation.  The  calculated 
momentum  thickness,  6,  shape  factor,  H,  boundary-layer  thickness,  6,  and 
shearing  stress,  t,  are  shown  in  Fig.  6.  B  ,  where  the  governing  equation  is 
changed  into  the  vorticity  equation,  is  chosen  at  2x/L  =  0.8  in  order  to  have 
reliable  values  for  initial  and  boundary  values. 

In  Figs.  7-9,  calculated  results  of  EM-125  are  shown. 


Fig.  7  shows  the  velocity  and  the  vorticity  distributions  at  2x/L  =  0.905 
and  1.10.  Calculated  results,  shown  at  every  0.1  L/U  time  step,  seem  to 
converge.  Agreement  with  the  measurements  is  not  so  good  at  2x/L  =  1.10  as 
that  at  2x/L  =  0.905. 

In  our  scheme,  we  imposed  the  symmetry  condition  on  y  =  ^ SYM ) *  We 
suspect  that  this  boundary  condition  may  be  the  cause  of  this  disagreement. 
As  is  seen  in  Fig.  9,  a  vortex  behind  the  body  is  predicted.  This  is  growing 
steadily  and  symmetrical ly.  Experimental ly ,  however,  it  is  observed  that 
vortices  are  unsteady  and  moving  across  the  symmetry  plane.  At  2x/L  =  0.905, 
flows  may  not  be  affected  so  much  by  the  symmetry  condition.  It  may  be 
necessary,  though  the  computing  time  and  the  storage  may  increase,  to  remove 
the  symmetry  condition  and  to  introduce  a  slight  disturbance  to  obtain  a  more 
realistic  flow. 

Fig.  8  shows  the  turbulence  quantities.  Undesirable  changes,  observed  in 
the  results  at  2x/L  =  0.905  around  2y/L  =  0.06,  are  after-effects  of  Cebeci- 
Smith's  model  which  is  invoked  to  determine  the  boundary  values  of  k 

and  e  on  B^. 

Fig.  9  shows  the  flow  patterns  at  three  time  steps.  Vectors  show  the 
velocity  and  the  flow  direction.  It  is  predicted  that  separation  is  occurring 
around  2x/L  =  0.975. 

Calculated  results  of  EM-200  are  shown  in  Figs.  10  and  11.  Fig.  10  shows 
the  computed  results  at  2x/L  =  0.95  and  1.10.  Here  the  calculations  were 
carried  out  up  to  t  =  0.4L/U.  A  much  longer  time  than  this  is  needed  for  the 
results  to  converge,  especially  far  downstream.  Indeed,  there  is  a 

possibility  that  the  results  may  not  converge.  Although  some  significant 
discrepancies  from  measurements  are  observed,  it  can  be  said  that  the  present 
method  predicts  the  near  wake  flow  fairly  well. 

Fig.  11  depicts  the  development  of  the  wake;  the  symmetric  vortex  becomes 
strong  and  eventually  a  significant  reverse  flow  is  realized.  The  separation 
position  is  predicted  to  occur  around  2x/L  =  0.925. 

These  calculations  for  the  elliptic  cylinders  raise  doubts  as  to  whether 

we  can  expect  a  steady  or  quasi-steady  flow  in  the  wake.  It  appears  to  be 

necessary  to  remove  the  symmetric  flow  condition  at  least. 
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5.  Concluding  remarks 

We  can  conclude  that  the  present  scheme,  where  the  time-dependent 
vorticity  transport  equation  is  used  in  its  full  form  together  with  the 
boundary- 1 ayer  equation,  is  promising.  Several  modifications  in  numerical 
techniques  seem  to  be  required;  e.g.,  the  prescribed  boundary  values  of  the 
vorticity  need  to  be  improved.  It  can  be  also  said  that  the  k-e  turbulence 
model  is  applicable  for  a  near  wake  prediction. 

Though  it  is  not  our  main  purpose  to  discuss  experiments,  it  is  highly 
necessary  to  carry  out  precise  velocity  measurements  in  the  near  wake  of  blunt 
bodies  which  can  provide  not  only  averaged  quantities  but  also  time-dependent 
quantities.  This  would  give  a  much  better  understanding  of  separated  flows 
and  accelerate  development  of  numerical  calculation  procedures. 
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Appendix 

Let  be  a  velocity  component  due  to  a  vorticity  distribution  in . 
Because  q  satisfies  the  continuity  equation,  it  can  be  expressed  in  terms  of 
a  vector  potential  E  in  the  form 


qy=  V  xE,  V.E  =  0 


(A— 1 ) 


A  curl  operation  applied  to  Eq.  (A-l)  then  yields 


!  E  =  -  «. 


(A-3) 


Thus  the  vector  function  E  is  a  solution  of  the  Poisson  equation  given  by  Eq. 
(A-3). 

By  making  use  of  Green's  theorem  in  the  domain  exterior  to  the  body,  we 
obtain 


4.E  ■///  ♦  //  ()|f  -  IS', 


(A-4) 


V  S< 


where  V  is  the  nonzero  vorticity  region  surrounded  by  S  and  S  which  are  a 

b  00 

body  surface  and  a  closing  surface  at  infinity  respectively;  n  denotes 
distance  along  the  outward  normal  to  sb  +  »  r  the  distance  between  the 

field  point  and  the  integrating  point;  primes  indicate  the  values  at  the 
latter  point. 

3 

If  ui  tends  to  zero  at  the  rate  of  1/R  ,  where  R  is  the  distance  to  , 
since  E  vanishes  as  1/R,  the  integration  over  becomes  infinitesimal ly  small 
when  R  is  chosen  large  enough.  Then  we  have 

4,1-111  fdV*//  (yff  -  E'  hf>dS'-  <A'5> 

M  SB 

Now  we  assume  another  vorticity  distribution  u>'  inside  Sg.  Similarly  to 
Eq.  (A-3)  we  have 


7  Ej=  -  a) j . 


(A-6) 


By  applying  Green's  theorem  to  the  interior  domain  of  Sn,  we  get 
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0-  111  or  -//  <}£■-  EjfjildS' 

B  Sg 

where  B  denotes  the  domain  interior  to  S_. 

D 

The  addition  of  Eq.  (A-7)  to  Eq.  (A-5)  yields 

I 

4itE  =  ///  f  dv*  +  ///  ~dV' 

V  B 

*//  (7k  <E'-Ei)  -  <e'-£T>!h7»  dS'- 

bB 

If  we  are  concerned  with  only  the  exterior  flow,  we  may  put 


(A-7) 


Er  E 


Yr  ?fr  -  Ei>  on  sb ’ 


(A-9) 


where  y^  is  a  circulation  distribution  whose  direction  is  tangential  to  S^. 
Then  Eq.  (A-8)  is  written  as 


«**-///  *  III  r  ♦  II  r  ds' 

V  8 


1 


(A— 10 ) 


Eq.  (A- 10 )  gives  a  solution  of  Eq.  (A-3).  Because  it  can  be  easily 
proved  that  E,  given  by  Eq.  (A-10),  satisfies  Eq.  (A-2),  the  induced  velocity 
is  given  by 

4n  qv=  III  (?  x  ~~ )  d*’  +  ///  (V  x  ^-)  d r  +JJ  (V  x  ~ )  dS1  (A- 1 1 ) 

V  B  SB 

Because  the  second  term  on  r.h.s.  of  Eq.  (A-ll)  has  a  potential  in  V  (this  can 
be  easily  proved  by  taking  the  curl  operation  to  it),  it  can  be  replaced  by  a 
source  distribution  over  SD.  Then  we  have  finally 

D 

a '  y  1 

///  (Vx^-)dT-vf  /  T^dS'  +  JJ  (Vx-pi-MS1,  (A-12) 

V  SB  SB 

where  Oj  and  are  determined  so  as  to  satisfy  the  hull  surface  condition. 

The  addition  of  the  velocity  component  given  by  gives  the  expression 
of  Eq.  (9)  in  the  text. 
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Table  2  Boundary  Conditions 
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Fig.l  Coordinate  System  and  Definitions 
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Fig.  10a  Comparisons  of  Velocity  Profiles  (EM-2^C.) 
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Fig.  10b  Comparisons  of  Velocity  Profiles  (EM- 
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C.  Free-Surface  Boundary  Layer  and  Necklace  Vortex  Formation 
1«  Introduction 

The  existence  of  a  shear  layer  beneath  a  wavy  free  surface  has  been 
experimental ly  suggested;  e.g.,  Takahama  [1]  and  Longuet-Higgins  [2], 
Batchelor  [3]  gives  a  theoretical  explanation  for  it c  development.  The  shear 
layer,  which  we  shall  call  a  free-surface  boundary  layer,  develops  due  to  the 
zero-shear-stress  condition  when  the  free  surface  has  a  significant  curvature. 

Besides  the  boundary  layer  on  the  hull  and  its  accompanying  wake,  a  free- 
surface  boundary  layer  may  also  affect  the  ship  wave  resistance.  The  viscous 
effect  of  the  free-surface  boundary  layer  on  wave  resistance  was  discussed  by 
Mori  [4],  wh.re,  by  order-of-magnitude  analysis,  it  was  concluded  that  its 
effect  might  be  less  than  those  of  the  boundary  layer  on  the  hull  and  the 
wake,  and  only  the  latter  was  included  in  his  calculation. 

Recently,  however.  Kayo  [5]  has  found  that  the  presence  of  a  free-surface 
shear  layer,  which  is  realized  artificially  in  his  research  signi ficantly 
affects  the  wave  resistance. 

In  the  first  half  of  the  present  report,  an  attempt  is  made  to  give 
equations  to  predict  the  free-surface  boundary  layer  in  the  2-D  case. 

In  the  latter  part,  a  theoretical  explanation  is  made  for  the  formation 
of  the  so-called  necklace  vortex  around  the  bow.  It  is  based  on  the  drift 
theory  (Lighthill  [6])  or  the  vortex-stretching  theory  (Batchelor  [7]).  A 
similar  idea  has  been  introduced  by  Baba  [8]. 

Through  a  simple  calculation  for  the  case  of  a  vertical  circular 
cylinder,  it  is  demonstrated  that  even  a  slight  vorticity  in  the  incident  flow 
can  become  very  intense  in  the  vicinity  of  the  bow  where  flows  are  retarded, 
and  that  the  observed  effect  of  a  bulbous  bow  on  inhibiting  the  formaiton  of 
this  vortex  can  be  explained  by  the  acceleration  of  the  flow  near  the  free 
surface  due  to  the  presence  of  the  bulb. 

The  combination  of  these  two  studies,  which  is  left  for  a  future 
investigation,  may  throw  light  on  the  bow-wave  phenomena. 
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2.  Free  surface  boundary  laurer  equations 

We  limit  ourselves  here  to  the  2-D  case.  An  orthogonal  curvilinear 
coordinate  system  (s,n)  is  used  together  with  the  Cartesian  coordinate  (o-xyz) 
system,  as  shown  in  Fig.  1.  The  normal  coordinate  n  consists  of  straight 
lines.  u  and  w  are  the  velocity  components  in  s-  and  n-directions 
respectively. 

The  momentum  equation  in  the  s-direction  is  given  by 


3  u 
h3  s 


3  U  . 

+  W  —  +  KUW 

3  n 


1  3  P 
p  h3s 


n  +  I  +  3t- 
y  h3s  p  'h3s  3n 


-  + 


2kt  ) , 


(1) 


where  h  is  the  metric  coefficient  of  the  s-coordinate,  P  the  pressure, 
p  the  density,  and  o  and  x  are  the  normal  stress  and  the  tangential  stress  in 
s-direction  respectively,  including  the  Reynolds  stresses.  k  is  the  curvature 
of  the  s-coordinate  defined  by 


< 


3h  _  ix 

han  .2 
dx 


{  1  + 


^dx;  ‘ 


-3/2 


(2) 


The  continuity  equation  is  given  by 


3  u 
h3s 


+ 


3w 

3n 


+  <w 


0. 


(3) 


We  assume  that  the  related  layer  is  so  thin  that  the  pressure  is  equal  to 
the  constant  atmospheric  pressure  and  that  the  derivative  with  respect  to  s  is 
smaller  than  that  with  respect  to  n.  Then  Eq.  (1)  can  be  approximated  by 


3  u 
h3s 


+  w 


3_u 

3n 


+  KUW 


=  -g  + 
y  h3  s 


1  /3x 
p  '3n 


2kt  ) , 


(4) 


where  c  is  the  free  surface  elevation  whose  gradient  is  written  in  the  steady 
case  as 


dc 

9  dt+  * 


43  . 

ds 


0, 


(5) 


where  q  is  the  potential  velocity. 

The  boundary  conditions  on  the  free  surface  and  at  the  edge  of  the 
boundary  layer  are  as  follows; 
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1 


n  =  0 


n  =  6 


t  =  0, 

u)  =  0,  u  =  q, 


where  w  is  the  vorticity  defined  by 


3  u  ,  3  w 

“  =  3n~  +  "u  '  hiss’* 


(6) 

(7) 


(8) 


We  assume  that  the  Reynolds  stresses  are  given  in  exactly  the  same  form 
as  the  stresses  due  to  the  molecular  viscosity  in  terms  of  velocity 
gradient.  Then  t  is  given  by 


T 


,3u 

=  “efo  '  KU  + 


3  w 
h3s 


), 


(9) 


where  Me  is  the  equivalent  viscosity  coefficient. 

If  the  flow  is  inviscid,  the  boundary  condition  of  Eq.  (6)  is 
automatical ly  satisfied  and  no  boundary  layers  develop.  If  the  flow  is 
viscous,  however,  Eq.  (6)  yields  an  additional  free-surface  condition 

§£  -  xu  =  0  ,  n  =  0  (10) 


which  is  the  source  of  a  possible  development  of  the  free-surface  boundary 
layer;  the  substitution  of  Eq.  (10)  into  Eq.  (8)  gives  the  vorticity  on  the 
free  surface; 


w  =  2icu  ,  n  -  0  (11) 

which  is  not  zero  unless  <  is  zero.  This  means  that  the  vorticity  is  not  zero 
on  a  curved  free  surface  and  a  boundary  layer  possibly  exists.  It  should  be 
remembered  that,  if  the  free  surface  is  flat,  the  boundary  condition  of  Eq. 
(10)  is  automatically  satisfied. 

In  the  present  paper,  an  integral  method  is  used  and  the  ordinary 
boundary-layer  approximations  are  invoked.  The  third  terms  on  l.h.s's  of  Eqs. 
(3)  and  (4)  are  omitted,  and  *  is  assumed  constant  across  the  boundary  layer 
and  equal  to  the  free-surface  curvature. 


| 


The  integration  of  Eq.  (4)  with  respect  to  n  yields 


ra?*  <2e 


+  6 


\  I  <*9_  = _ L. 

'  q  hds  2 

H  pq 


l  n=6 


2<_ 


J  xdn. 


(12) 


where  9  and  <$*  are  the  momentum  thickness  and  the  displacement  thickness 
respectively,  defined  by 


q29 


5 

J  u(q-u)dn  , 
o 


6 

q6  *  =  j  (q-u)dn. 
o 


(13) 


In  the  derivation  of  Eq.  (12),  we  have  used  Eqs.  (3)  and  (5). 

After  neglecting  the  third  term  of  Eq.  (9),  its  substitution  into  Eq. 
(12)  yields 


d8 

hds 


-  (ze  *  <*>!&*  U  <0 


q  hds 


2%» 


(14) 


where  qQ  is  the  velocity  on  the  free  surface. 
If  Eq.  (8)  can  be  approximated  as 


3u  A  n 

10  =  Fn  +  KU  =  0 


(35) 


at  n  =  6  ,  the  r.h.s.  of  Eq.  (14)  can  be  simplified  to 


<z»  *  «*>  5  sli  ‘  z* e' 

As  a  second  equation,  we  use  the  entrainment  equation*) 
by  an  integration  of  Eq.  (3), 


where 


d 

hds 


(HEe)  +  HEe 


I  <3a._  _  F 

q  hds  L 


(16) 

which  is  obtained 

(17) 

(18) 


*)of  course,  we  can  use  an  alternative  equation  such  as  the  moment  of  moment 
equation.  The  choice  of  the  best  one  is  left  for  future  work. 
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The  entrainment  function  E  is  given  apriori  as  a  function  of  H^. 

Now  we  express  the  velocity  profile  in  the  boundary  layer  by  a  polynomial 
summation 


u/q  =  a0+  ai(£)  +  a2(J)2  +  a,(j)3 


^6 


(20) 


Then  the  boundary  conditions  Eqs.  (6)  and  (7)  give 


al  =  kV 


d2  ”  -a0(2k  +  3)  +  k  +  3> 


(21) 


a3  =  (aQ-l)(k  +  2), 

where  k  =  k6  ;  the  approximation  of  Eq.  (15)  is  also  invoked. 

Then  &*  and  e  are  written  in  terms  of  a  and  k  as 

o 


6*/6  =  '  12  {ao  (k  +  6)  +  (k"6)}  ’ 

6/6  =  -  {  8a2(k2  +  Ilk  +  39)  +  2aQ(6k2  +  17k  -  102) 


(22) 


+  8k  +  18k  -  108)}. 


Once  e  and  have  been  obtained  at  a  certain  position  by  solving  Eqs. 
(16)  and  (17),  then  aQ  and  k  can  be  determined  from  Eqs.  (18)  and  (22). 

Thus  we  can  calculate  the  free-surface  boundary-layer  flow.  Further  we 
can  include  it  as  one  of  the  double-hull  velocity  components  in  the  Rankine 
source  method  (see  part  A). 


Formation  of  necklace  vortex 

We  assume  that  the  neck  lac  vortex  around  the  bow  is  formed  by  an 
accumulation  of  the  vorticity  generated  in  the  free-surface  boundary  layer. 
As  the  drift  theory  suggests  [6],  because  of  the  stretching  of  the  streamlines 
around  the  bow,  even  a  slight  vorticity  can  form  such  an  intensive  vortex  as 
is  observed. 
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In  the  inviscid  fluid,  a  vortex  tube  moves  with  a  fluid  particle  and  its 
strength  remains  constant  (Helmholtz's  theorem).  This  theorem  gives  a 
relation 


li) 


It+Atl  =  &lt  +  At) 

u>(trf  juTtir 


(23) 


where  i  is  a  material  element  parallel  to  the  vorticity  u>  [7],  and  t  is  time. 

To  solve  Eq.  (23)  is  equivalent  to  solving  the  vorticity  transport 
equation  with  the  viscous  diffusion  term  neglected. 

If  the  vorticity  can  be  assumed  to  be  distributed  on  the  z  =  0  plane,  the 
induced  velocity  is  given  by 


1  . 

V  2  “x5’ 

\  --K«  ■  (24> 

«,  ■  -  L  II  S  i*  -  *')»v-  (y  -  y')  «•„>  d*'  dy'- 

r  J 

where  r2  =  (x  -  x' )2  +  (y  -  y')2,  u>v  and  w  are  x- ,  and  y-  components 

x  y 

of  L>  respectively.  6  is  the  nonzero  vorticity  layer.  Here  we  assume  that 
the  vorticity  is  constant  in  the  z-direction  within  the  layer. 

The  velocity  given  by  Eq.  (24)  can  be  included  in  qy  in  the  Rankine 
source  method  (see  Part  A),  and  we  can  calculate  the  free-surface  elevation 
including  the  effect  of  a  necklace  vortex. 

By  making  use  of  the  relation 

gVdH  =  U  X  q,  (25) 


we  obtain  the  headloss  6H  as 


x 

g6H  =  J  ui  wdx 

x  y 
o 

In  Eq.  (26),  xq  is  a  far-upstream  position  where  the  head  loss  is  zero. 


(26) 
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L  A  Nugerical  example  for  a  vertical  cylinder 

When  a  vertical  cylinder  is  in  a  uniform  flow  U,  the  velocity  field  is 
given  by 


u  =  U 


v 


(27) 


where  a  is  the  radius  of  the  cylinder. 

Fig.  2  shows  the  free-surface  elevation  along  the  centerplane  (y  =  0) 
which  is  approximately  obtained  by  making  use  of  Bernoulli's  formula  (the 
kinematic  condition  is  not  satisfied).  The  measurement  was  carried  out  with  a 
point  gauge.  Agreement  is  good  except  near  the  bow.  Just  in  front  of  the 
bow,  a  vortical  motion  can  be  observed.  Around  x/a  =  -2.1,  a  sudden  increment 
in  the  surface  elevation,  which  seems  to  form  a  wave-front  line,  is  seen. 

Fig.  3  shows  the  paths  of  fluid  material  particles,  obtained  by 

t 

x  +  iy  =  xQ  +  iyQ  +  j  (u  +  iv)  dt,  (28) 

o 

where  (xQ,yo)  is  the  initial  position  of  a  particle  along  which  a  vortex  tube 
is  initially  assumed.  In  our  present  calculation,  it  is  chosen 
at  x/a  =  -1.9.  The  initial  value  of  vorticity  is  determined  by  Eq.  (11) 
as  |w  |  =  1.0. 

The  solid  lines  show  streamlines,  the  broken  lines  the  positions  of 
particles  at  the  related  time.  It  can  be  observed  that  vortex  tubes  are 

greatly  stretched  in  the  vicinity  of  the  cylinder. 

Fig.  4  shows  the  equi vorticity  contours;  Fig.  4a  and  Fig.  4b  show  those 

of  w  and  a)  respectively.  In  front  of  the  cylinder,  <u  has  increased  more 

x  y  y 

than  seven-times  its  original  intensity.  This  accumulation  may  explain  the 
existence  of  the  experimentally  observed  vortical  motion  there. 

The  effect  on  the  generation  of  the  component  is  drastic.  Although 

the  <*>x-  component  rarely  exists  along  the  initial  distribution  line,  a 

significant  vorticity  is  attained  which  may  form  the  longitudinal  vortex  along 
the  cylinder. 


Such  accumulations  may  possibly  affect  the  velocity  field.  Fig.  5  shows 
the  induced  velocity  due  to  the  free-surface  vorticity  which  is  obtained  from 
Eq.  (24).  6  is  assigned  a  constant  value  of  0.02a  over  the  whole  domain  which 
is  roughly  equal  to  10%  of  the  wave  height.  Although  the  values  assigned  in 
the  present  calculation  are  quite  modest,  the  induced  velocity  is 
significant.  The  negative  sign  of  uv  indicates  that  the  induced  velocity  is 
backward.  As  seen  in  Fig.  5,  a  stagnation  point,  where  the  addition  of  uv  to 
u  becomes  zero,  is  realized  around  x/a  =  -  1.05.  This  may  account  for  the 
observed  vortical  motion  ahead  of  the  bow,  and  possibly  for  the  occurrence  of 
wave  breaking.  An  exact  calculation  of  the  free-surface  elevation,  including 
the  induced  velocity,  is  highly  desirable. 

In  Fig.  5,  the  head  loss,  obtained  from  Eq.  (26),  is  also  shown. 
Although  it  does  not  always  compensate  all  the  discrepancies  observed  in  Fig. 
2,  it  may  be  enough  to  produce  an  increment  in  resistance. 

It  is  well  known  that  a  bulbous  bow  reduces  resistance  even  at  a  low 
Froude  number.  The  reason  is  believed  to  be  that  bulbs  improve  the  bow  flow 
so  as  to  inhibit  secondary  flows  or  vortex  formation  but  this  has  not  been 
clearly  demonstrated.  A  simple  attempt  will  now  be  made  to  shed  light  on  this 
phenomenon  by  introducing  a  point  doublet  in  front  of  the  cylinder.  A  point 
doublet  of  moment  strength  0.125  Ua3  is  submerged  at  (-1.5a,  0,  -0.75a).  The 
original  flow  field,  given  by  Eq.  (27),  does  not  correspond  to  that  of  the 
circular  cylinder  any  more;  this  is  ignored. 

Fig,  6  shows  the  paths  of  material  particles.  Due  to  the  acceleration  by 
the  point  doublet,  the  stretching  of  a  material  tube  is  greatly  moderated. 
Fig.  7  shows  the  equi vorticity  contours;  as  observed,  the  accumulation  of 
vorticity  is  much  less  than  that  without  a  doublet.  The  longitudinal  vortex 
has  been  remarkably  weakened.  The  induced  velocity,  shown  in  Fig.  5,  is  a’so 
much  less. 

These  results  show  that  bow  bulbs  play  a  role  to  make  the  accumulation  of 
vorticity  less.  This  fact  may  provide  a  possible  guideline  to  hull-form 
design;  bulbs  should  be  designed  in  order  for  the  free-surface  flow  (not  the 
entire  flow!)  to  be  accelerated  so  as  to  avoid  the  intensification  of  the 
free-surface  vorticity. 


4.  Concluding  remarks 

The  boundary  layer  equations  are  derived  and  it  is  shown  that  a  free- 
surface  boundary  layer  can  develop  when  a  free  surface  has  a  significant 
curvature. 

The  calculations  and  experiments,  which  will  provide  basic  data  for  this 
not-well -known  problem,  are  being  performed. 

The  simple  theory  of  drift  is  applied  to  the  free-surface  vorticity.  The 
formations  of  the  necklace  vortex  and  the  longitudinal  vortex  along  the  hull 
have  been  demonstrated.  A  thorough  computation  to  predict  the  free-surface 
elevation  will  be  undertaken  by  combining  the  described  procedures  with 
Rankine  source  method. 
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Fig.  7a  Equi-Vorticity  Contour  around  Cylinder 

with  Submerged  Point  Doublet  (^x) 


